Effective spin model for the spin-liquid phase of the Hubbard model on the triangular lattice 
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We show that the spin liquid phase of the half-tilled Hubbard model on the triangular lattice can be described 
by a pure spin model. This is based on a high-order strong coupling expansion (up to order 12) using perturbative 
continuous unitary transformations. The resulting spin model is consistent with a transition from three-sublattice 
long-range magnetic order to an insulating spin liquid phase, and with a jump of the double occupancy at the 
transition. Exact diagonalizations of both models show that the effective spin model is quantitatively accurate 
well into the spin liquid phase, and a comparison with the Gutzwiller projected Fermi sea suggests a gapless 
spectrum and a spinon Fermi surface. 
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Although the Hubbard model has been one of the central 
paradigms in the field of strongly correlated systems for about 
five decades, new aspects of its extremely rich phase diagram 
are regularly unveiled. Even at half-filling, the popular wis- 
dom according to which the model has only two phases, a 
metallic one at weak coupling and an insulating one at strong 
coupling separated by a first order transition [1], has been re- 
cently challenged. This goes back to the work of Morita et 
at. [2] on the triangular lattice which revealed the presence of 
a non-magnetic insulating phase close to the metal-insulator 
transition using path integral renormalization group. Further 
evidence in favour of a phase transition inside the insulating 
phase has been reported using a variety of theoretical tools [3- 
6]. More recently, an intermediate spin liquid (SL) phase has 
also been identified on the honeycomb lattice using Quantum 
Monte Carlo simulations [7]. 

The precise nature of the SL phase of the Hubbard model 
on the triangular lattice is of direct experimental relevance for 
the 2D organic salt k-(BEDT-TTF)2Cu2(CN)3 [8]. As such, it 
has already attracted a lot of attention, but fundamental ques- 
tions such as the appropriate low-energy effective theory re- 
main unanswered. Since the phase is insulating, an effective 
model where charge fluctuations are treated as virtual excita- 
tions should be possible. One step forward in this direction 
has been taken by Motrunich [9], who proposed to describe 
the SL phase with 4-spin interactions. However, whether a de- 
scription in terms of a pure spin model is possible is far from 
obvious, in particular since there seems to be a jump in the 
double occupancy at the transition from the three-sublattice 
Neel phase to the SL [2, 6]. 

In this Letter, we show that the correct low-energy theory of 
both insulating phases, and in particular of the SL phase, is in- 
deed a pure spin model. This has been achieved by deriving an 
effective spin model to very high order about the strong cou- 
pling limit using perturbative continuous unitary transforma- 
tions (PCUTs), and by showing that it gives a qualitative and 
quantitative account of the transition from the three-sublattice 
magnetic order to the SL state. This description gives deep 
insight into the nature of the SL phase and clearly provides 
the appropriate framework for further studies. 



The starting point is the single-band Hubbard model on the 
triangular lattice defined by the Hamiltonian 



n = Hu + Ht 
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where the sum over i runs over the sites of a triangular lat- 
tice, and the sum over over pairs of nearest neigh- 
bors. To derive an effective low-energy Hamiltonian, we use 
PCUTs [10-13] about the localized limit treating the hopping 
term as a perturbation. The kinetic part can be written as 
Ht = t (To + Ti + where Tm changes the number of 
doubly-occupied sites by m. The PCUT method provides or- 
der by order in t/U an effective Hamiltonian Hcff with the 
property [Hu,Hcs] = 0, i.e. Hcff is block-diagonal in the 
number of doubly-occupied sites. The low-energy physics at 
half filling is expected to be contained in the block with no 
doubly-occupied site H^ff, which can be expressed as a spin 
1/2 model of the form 

H^s = const + J2jfi {^r ■ Sf+n) (2) 

r,n 

r,ni ,n2 ,^3 

where r and f+fia denote sites on the triangular lattice. All 
remaining terms can be written in a similar way as products 
of S,f ■ Sf due to the SU(2) symmetry of the Hubbard model. 

The PCUT provides the magnetic exchange couplings as 
series expansions in t/U. Since the spectrum of the Hub- 
bard model is symmetric at half filling under the exchange 
t o —t, only even order contributions are present [14]. We 
have determined all 2-spin, 4-spin and 6-spin interactions up 
to order 12. The obtained series are valid in the thermody- 
namic limit. To this end we have fully exploited for the first 
time in a PCUT calculation the linked cluster theorem by us- 
ing graph theory. At order 12 there are 1336 topologically 
different graphs. The major complication for the determina- 
tion of higher orders comes from the large number of different 
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FIG. 1: (Color online) The most relevant exchange couplings dis- 
played as a function of t/U. The 12th order bare series of all cou- 
plings appearing up to order 4 and the largest 4-spin interaction Jq^ 
plus the largest 6-spin interaction J^'^^ appearing in order 6 are plot- 
ted. Grey shaded area marks the region where extrapolation effects 
become important. Vertical dashed line indicates the magnetic phase 
transition. Inset: various Pade approximants as a function of t/U for 
the nearast-neighbor Heisenberg exchange Ji_o. 



FIG. 2: (Color online) (a) Energy per site as a function of U /t for the 
Hubbard model and the effective spin model, (b) Level crossings in 
the effective spin model signalled by the energy difference of the two 
lowest energies, (c) Double occupancy (times L*'^) as a function of 
U/t. (d) Magnetic structure factor at the 120° AFM ordering wave 
vector. Dashed lines denotes results for a simplified 2B-4B model 
evaluated up to order ten. Continuous lines with symbols are the 
results for the full 2B-4B-6B model considered here. The observable 
is up to order two in both cases. 



spin operators, which requires to calculate all possible matrix 
elements for each graph. 

The most relevant spin couplings are shown in Fig. 1 as a 
function of t/U. As already pointed out in Ref. 9, the most 
important correction to the nearest-neighbor Heisenberg ex- 
change at order 4 are the 4-spin interactions on a diamond 
plaquette. High-order contributions result in more complex 
multi-spin interactions with larger spatial extension but also 
lead to renormalizations of already existing couplings. The 
latter effect causes for example a sizable splitting of the two 
4-spin interactions shown in Fig. 1, which have the same ab- 
solute value at order 4. Processes involving sites far apart 
have small amplitudes, and multi-spin interactions between 
sites not far from each other are clearly the dominant correc- 
tions to the spin operators at order 4. 

At large U/t, nearest-neighbour exchange dominates, and 
the system is expected to develop three-sublattice long-range 
order [15-17]. Using exact diagonalizations (ED) of clusters 
up to 36 sites, we have investigated how the properties of the 
model evolve upon reducing U/t. The most striking feature 
is a level crossing in the ground state (GS) at U/t ~ 10 (see 
Fig. 2). This level crossing is observed for all cluster sizes at 
essentially the same value, which indicates the presence of a 
first order transition in the thermodynamic limit. Since this 
ratio is comparable to the critical ratio where other methods 
have detected a transition into a SL phase, it seems plausible 
that this level crossing corresponds to this transition. Further 
investigations below U/t ~ 10 confirm this guess (see be- 
low), but let us first investigate to which extent the effective 
spin model provides an accurate description in this parameter 
range. The proof that this is so relies on four observations: 



- Convergence of the series: As usual, we use Pade ap- 
proximants to extrapolate the series, and we consider that 
series are well converged as long as different Pade approxi- 
mants give consistent results. Here this is definitely true up to 
t/U ~ 0.125 (U/t ~ 8), i.e. well beyond the level crossing 
(see inset of Fig. 1). For all couplings, the absolute change 
between order 10 and 12 is below 10""* (10"^) foj. ^/jj = q.I 
{t/U ~ 0.125). In fact, clear indications of divergence only 
appear above t/U ~ 0.15, and this is probably related to the 
metal-insulator transition. 

- Comparison of GS energy with Hubbard model: In view 
of the anomaly observed at the transition by other approaches 
[2, 6], it is legitimate to ask whether the GS remains in the 
sector with no doubly-occupied site. To address this point, 
we have used ED of finite clusters to calculate the GS en- 
ergy of both the Hubbard model and the effective spin model 
as a function of U/t (see Figs. 2a and 2b) up to 21 sites for 
the Hubbard model, and up to 36 sites for the effective spin 
model. For the spin model, the total number of spin opera- 
tors is enormous, and truncating the hierarchy of operators is 
necessary. We have included in our calculation all spin oper- 
ators appearing up to order 6. Additionally, we have omitted 
for a given ratio t/U all terms whose coupling was smaller 
than 10^^. When embedded on our largest 36 site cluster, this 
corresponds to ^ 16'000 distinct terms in total. The resulting 
total energies per site are displayed in Fig. 2. The energy per 
site of the effective spin model depends very little on the size, 
and for the 21 -site cluster, the energies of the Hubbard model 
and of the effective spin model are in very good agreement on 
both sides of the transition. This last observation clearly es- 
tablishes that the effective model remains accurate in the SL 
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FIG. 3: (Color online) Comparison of the low-energy spectrum of 
a A'^s =21 Hubbard model al U/t — 9 (left panel) to results ob- 
tained with the effective spin model (right panel). The a::-axis labels 
the different symmetry sectors, while empty (filled) symbols denote 
= 1/2 (S^ — 3/2) levels, a, b and c label momenta in the inte- 
rior of the Brillouin zone. 



FIG. 4: (Color online) ED spectra of the truncated spin model for (a) 
Ns = 21, (b) iV^ = 27, and (c) iV^ = 36 as a function of U/t. Open 
(black) symbols denote non-magnetic levels. Filled (red) symbols 
are energy levels corresponding to excitations with finite spin above 
the GS. The grey shaded area signals the extension of the SL phase. 



phase. The GS energies for both the Hubbard model and the 
effective model are actually significantly lower than recent es- 
timates based on PIRG [6], see Fig. 2a. 

- Double occupancy: The comparison of the GS energies 
suggests that the GS is in the sector of the effective model 
with no doubly-occupied site in the U/t range considered. 
This result seems to be incompatible with the jump of the 
double occupancy Ud reported by other methods for the Hub- 
bard model [2, 6], but in fact it is not because, when calcu- 
lating the expectation value of observables from the effec- 
tive model, one must apply the same unitary transformation 
[12, 18] which, in the present case, leads to a non-zero value 
of Ud- In practice, it is simpler to use the Feynman-Hellmann 
theorem which allows to deduce it directly from the energy per 
site as Ud = du{TL)/N. The results are depicted in Fig. 2c. 
A jump is clearly present at the transition for all finite-size 
samples as a consequence of the level crossing of the effective 
model, and its magnitude is consistent with that calculated for 
the Hubbard model on the 21-site cluster (see Fig. 2c). The 
significant difference between 21 and 27 sites does not allow 
to perform a reliable finite-size scaling, but this does not af- 
fect the conclusion that the effective spin model captures the 
double occupancy jump at the transition. 

- Low-energy spectrum and GS momentum: It is also im- 
portant to check how well the effective model reproduces the 
low-energy spectrum of the original model. To address this 
issue, we compare the low-energy spectrum at U/t = 9 of 
a Ns = 21 Hubbard model with the low-energy excitations 
of the spin model on the same cluster [25] in Fig. 3. Quali- 
tatively, the agreement is good. In particular, the non-trivial 
symmetry sectors of the GS and of the first excited states are 
correctly reproduced. Quantitatively, the difference between 
the energies (a fraction of a percent) is dominated by the fact 
that the effective model on a finite cluster should include, 
starting from order 6, a few processes which wrap around 
the cluster These processes are not taken into account in our 



effective model since they are not present in the thermody- 
namic limit. The systematic error in the thermodynamic limit 
is much smaller since it only comes from the precision of the 
coefficients and the neglect of truncated terms. 

The truncated spin model is thus an excellent starting point 
to investigate the quantum magnetism of the Hubbard model 
on the triangular lattice. Consequently, we turn to a detailed 
discussion of its physical properties as a function of U/t. The 
low-energy ED spectra for clusters of Ng = 21, 27, and 36 
sites are shown in Fig. 4. For clarity the lowest level at each 
U/t value is taken as the reference energy and set to zero, and 
the spatial symmetry quantum numbers of the levels are not 
specified. For the 36-site cluster, the calculation of the excita- 
tion spectrum is numerically very expensive, and it has been 
performed only for three values of U /t close to the transition, 
and for a simplified two- and four-spin coupling model (in 
contrast to the ground state energy calculations in Fig. 2(a) 
which have been performed based on the more demanding 
model with the aforementioned cutoff 10^^). At large U/t 
(down to U/t ~ 15), the level structure revealed by the 21 
and 27-site clusters is typical of the tower of states of a 120° 
antiferromagnet (AFM) [15]. For smaller U/t ratios, some 
levels start to bend down and to deviate from the leading 1/U 
behavior of the Heisenberg regime. Drastic changes occur in 
the spectrum at U/t ?s 10. First of all, on all clusters, ex- 
cited states with different quantum numbers cross the large 
U/ 1 GS level. Besides, a large number of singlets accumulate 
at low energy in both the effective spin model and the Hub- 
bard model for U/t < 10. 

To better characterize this change of behaviour, we have 
calculated the structure factor corresponding to the three- 
sublattice 120° order as a function of t/U by implementing 
the unitary transformation on the observable at order t^ /U"^. 
Above the level crossing, this structure factor increases with 
the size, as expected for a long-range ordered phase. By con- 
trast, below the level crossing - i.e. in the putative SL region 
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- the structure factor is finite size saturated, implying either 
short range correlations, or an algebraic decay with a size- 
able decay exponent. At large U /t, the slight reduction of the 
structure factor is mainly a consequence of the admixture of 
empty and doubly occupied sites, which reduces the effective 
local spin length below 3/4. However, as we approach the 
level crossing, the extra terms in the spin Hamiltonian further 
reduce the structure factor, as can be seen by comparing the 
results for the 'full' effective model with those for a truncated 
model where only the 2-spin and 4-spin terms that appear at 
fourth order (but evaluated using the full series) are kept (see 
Fig. 2d). So one cannot exclude that the long-range magnetic 
order disappears slightly before the level crossing. 

The issue of the spin gap in the SL phase requires a careful 
discussion. Indeed, the energy of the first magnetic excitation 
is nearly constant between 21 and 27 sites and drops dramat- 
ically for 36 sites. So depending upon whether one includes 
the 36 cluster or not, finite-size effects suggest a vanishing 
gap or a finite gap. Now, the significant difference between 
the 36-site cluster and the other two might be due to an even- 
odd Ns effect, or to the use of a simplified model, and it might 
be better not to include it in the finite-size analysis. Now, con- 
centrating on the 21 and 27 site clusters, the energy of the 
first magnetic excitation is indeed almost constant, but if one 
looks carefully at the drastic U /t dependence of the excited 
magnetic levels when approaching the SL regime from large 
U /t, some excitations bend down, and their energy decreases 
very significantly between 21 and 27 sites, possibly leading 
to a gap closing in the thermodynamic limit. So it is difficult 
to decide on the presence of a spin gap based on the available 
ED data alone. In the next paragraph we show that the data 
can be interpreted from a broader perspective. 

In his seminal paper Motrunich [9] proposed a plain 
Gutzwiller projected Fermi sea as a prototype spin wave func- 
tion for the SL phase on the triangular lattice. This idea 
later materialized into a theory of Spin Bose Metal (SBM) 
phases [20] (see Ref. [21] for a related scenario), which 
has been successfully applied to a one-dimensional triangu- 
lar strip model with two- and four-spin interactions [20, 22]. 
Interestingly, the numerical low-energy spectrum of the trian- 
gular strip in the SBM phase looks qualitatively similar to the 
spectrum we observe on the triangular lattice, with respect to 
both the dense excitation spectrum with many low energy sin- 
glets and the irregular finite size behavior of the spin gap [22]. 
A second remarkable agreement between the SBM picture by 
Motrunich and our ED analysis comes from the quantum num- 
bers of the ground state wave function in the SL region. On 
the Ns = 27 site sample for example the noninteracting half 
filled Fermi sea has a sixfold degenerate ground state momen- 
tum. This prediction matches precisely the ground state quan- 
tum number of the effective spin model in the SL region adja- 
cent to the magnetically ordered phase. On the other samples 
{Ng = 21, 36) the Fermi sea ground state degeneracy is larger 
and thus less constraining, but the ground state quantum num- 
bers of the SL regime are still compatible with this picture. 

In conclusion, we have derived the appropriate low-energy 



theory of the non-magnetic phase of the half-filled Hubbard 
model on the triangular lattice. It is a pure spin model with 
multi-spin exchange. Furthermore, we have performed an ED 
investigation of this model. Several aspects of the results are 
compatible with the SBM picture predicting a gapless spec- 
trum with a spinon Fermi surface. It will be interesting to 
extend our investigation to the anisotropic triangular lattice 
[23, 24], and to compare the results directly to the SL mate- 
rial K-(BEDT-TTF)2Cu2(CN)3 and to the other members of 
the family. Finally, this approach is also applicable to other, 
possibly frustrated lattice topologies for which similarly inter- 
esting Mott phases can be expected. 
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